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Abstract. The coupled in-plane diffusion dynamics between point-particles 
embedded in stacked fluid membranes are investigated. We calculate the contributions 
to the coupling longitudinal and transverse diffusion coefficients due to particle motion 
within the different as well as the same membranes. The stacked geometry leads to a 
hydrodynamic coupling between the two membranes. 

Biological membranes are fundamental to the existence of life with their ability to 
separate the "in" and "out" of a cell. Typical membranes are composed of lipid molecules 
which can spontaneously self-assemble into a fluid bilayer structure pQ. Under normal 
physiological conditions, finite temperature induces Brownian motion of membrane 
constituents, resulting in diffusive transport. In the seminal paper by Saffman and 
Delbriick [2], the diffusion coefficient of a rigid protein in a membrane was calculated. 
Following this work, there have been theoretical studies on diffusion of a rigid disk on 
supported membranes [3] , a liquid domain in monolayers [U [5] , and a rod on immersed 
membranes [6j [7] or on Langmuir monolayers [8] . 

The present work focuses on the coupling of diffusion dynamics between two fluid 
membranes through an intervening bulk fluid. It was reported that model experimental 
systems in which two lipid membranes are stacked on a substrate exhibit correlated 
dynamics [9]. This planar geometry has become favorable to study the membrane 
dynamics which are otherwise not possible in vesicles [10]. The coupling effect between 
two membranes can be important in biological systems with large concentration of 
cells such as in tissues. Other examples are Gram-negative bacteria which enclose 
a periplasmic space with an approximate width of about 15 to 20 nm between their 
inner and outer lipid bilayers [TJ. Highly folded membranous organelles such as Golgi 
apparatus also correspond to a situation in which membranes come in close proximity to 
each other. In all these cases, it is very relevant to consider the hydrodynamic coupling 
between two biomembranes. 

Hydrodynamic models of membranes exploit their fluid nature, in which the 
membranes are assumed to be two-dimensional (2D) viscous fluid sheets embedded in a 
three-dimensional (3D) solvent. One such investigation involved the calculation of the 
correlated diffusion of proteins embedded in a membrane immersed in an unbounded 
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Figure 1. Schematic picture showing a set of stacked fluid membranes (solid thick 
lines) having 2D viscosities r\ located at z = ±h/2. The two membranes are denoted 
by the labels 1 and 2. The solvent of 3D viscosities r] s are labeled by the regions i, ii 
and hi. 



fluid [TT] or a membrane adjacent to a solid support [12]. In this Communication, we 
report on the coupling diffusion coefficients of particles embedded in such a stacked 
membrane system. 

We first establish the governing equations. It is assumed here that the membranes 
are infinite planar sheets of fluid. The persistence length is defined as the length over 
which the correlations of the normal vector decay to zero [13]. For a lipid bilayer 
membrane with bending rigidity of about 10 k^T, the persistence length turns out to be 
much larger than the size of a typical cell. Hence we neglect out-of-plane fluctuations 
of the membrane for mathematical simplicity. The fluid membranes are embedded in 
a bulk solvent such as water or suitable buffer solution. As shown in figure [IJ the 
membranes are fixed in the xy-plane at z = ±h/2. Let v^r) be the 2D velocity of 
the membrane fluids. Here the index i — 1,2 represents the two membranes, and the 
2D vector r = (x, y) represents a point in the planes of the membranes. We work in 
the low-Reynolds number regime of the membrane hydrodynamics so that the inertial 
effects can be neglected. This allows us to use the 2D Stokes equations 

V V 2 V® - + f s (4) + F w = 0, (1) 

along with the incompressibility condition, 

V • v (i) = 0. (2) 

Here V is a 2D differential operator, rj is the 2D membrane viscosity (same for both the 
membranes), p^(r) the 2D in-plane pressure, f s (r) the force exerted on the membrane 
by the surrounding fluid, and F^^r) is any other force acting on the membrane. 
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The solvent regions are denoted by the index j = i, ii, iii. The velocities and 
pressures in these regions are written as v^(r, z) and p^>(r, z), respectively. We assume 
that the solvent in the three regions have the same 3D viscosity denoted by r] s . The 
solvent inertia is neglected and hence it also obeys the 3D Stokes equations 

Vb V 2 v U) _ Vj9 (i) = 0, (3) 

where V represents a 3D differential operator. Similar to the fluid membrane, the solvent 
in all the regions are considered to be incompressible 

V ■ v (i) =0. (4) 

The presence of the surrounding solvent is important because it exerts force on the 
liquid membranes. The force on membrane 1, indicated as f s in (jTJ), is given by the 
projection of (cr^ — <r^) z= h/2 • &z onto the xy-p\ane of the membrane. Here e z is the 
unit vector along the z-axis, and cr^ are the stress tensors 



Vv« + (Vv«) T . (5) 



In the above, I is the identity tensor and the superscript "T" indicates the transpose. 
Similarly, the force on membrane 2, is given by the projection of (cr^ — (T^) z =-h/2 • e z 
on the xy-plane. The general procedure is to first resolve ([1]) into components along k 
and perpendicular to it, where k is the 2D wave vector in Fourier space. We then solve 
the resulting differential equations for the velocities. Stick boundary conditions at the 
membrane-solvent interfaces are imposed. It is also assumed that the solvent velocities 
decay to zero at sufficiently large distances from the membranes. A similar calculation 
for a single membrane confined between two walls has been previously performed [TJJ[T5]. 

Owing to the linearity of governing Stokes equations, the in-plane velocity in 
membrane 1 can be obtained in Fourier space as 

^ M = G<£> [k] if> [k] + G^f [k] Ff [k] . (6) 
Here G^ and G^ (a, (3 = x,y) are the mobility tensors given by 
,(11),., 1 1 + 2K + coth(KH) ( kjza 



°^ [k] "^ Kg(K,H) (**-iH' (7) 



( i2) 1 cosech(KH) f k a k p » 

°^ [k] = ^ Kg(K,H) l^- — )' (8) 



with 



g{K, H) = l + 2K(1 + K) + (1 + 2K) coth(M), (9) 

and v = 2r] s /r]. Notice that v~ l represents the hydrodynamic screening length. We have 
also used the definitions H = hv and K = kjv with k = |k|. By symmetry, a similar 
set of expressions can be written down for the membrane 2 also. The above mobility 
tensors represent the hydrodynamic coupling between the two membranes. 

Consider a pair of particles embedded in the membrane undergoing Brownian 
motion separated by r. For sufficiently large enough times, the particle displacements 
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Figure 2. Dimensionless Df^ (R,H) as a function of R for various values of H. 
The red circles, blue squares and green triangles correspond to H = 0.1, 1 and 
10, respectively (color online). The solid line corresponds to (jlll) representing the 
analytical expression for large H limit. The dashed line corresponds to the small H 
limit. 



obey (Ar a Ar^) = 2D a p{r)t where Ar a represents the displacement of the first particle 
and Ar'p represents that of the second particle. In the above, the diffusion tensor for 
over-damped dynamics is given by the Einstein relation D a p = k^TGap where k B is 
the Boltzmann constant and T the temperature. We now apply these definitions to 
a double-membrane system. The line of centers connecting any two particles in the 
membranes after projection on to the 2D plane can be taken to be along the x-axis 
without loss of generality. Then we obtain the coupling longitudinal diffusion coefficients 
as (re x ) of two particles within the same membrane, and (r) = 

ksTGx^irex) of two particles in different membranes. The coupling transverse 
diffusion coefficients are D^\r) = k B TG { yy\re x ) and D$ 2 \r) = k^Gyf (reA The 
longitudinal coupling diffusion coefficient is associated with Brownian motion along the 
line of centers, while the transverse one is associated with motion perpendicular to 
the line of centers [TT1 [T2| IT6] . In the following, we discuss the above four diffusion 
coefficients sequentially. 

The real-space expressions for the mobility tensors can be obtained by inverse 
Fourier transform. Using the notation R = ru, we obtain the longitudinal coupling 
diffusion coefficient within the same membrane as 



where J n (t) are the Bessel function of the first kind of order n. Hereafter, all the diffusion 
coefficients are scaled by a factor k B T / (Anri) in order to make them dimensionless. In 
figure [21 we plot d£ 1] (R,H) CIS db function of R for various H following numerical 




(10) 
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Figure 3. Dimensionless D^\ll,H) as a function of R for various values of H. 
The red circles, blue squares and green triangles correspond to H = 0.1, 1 and 
10, respectively (color online). The solid line corresponds to (|13l) representing the 
analytical expression for large H limit. The dashed line corresponds to the small H 
limit. 



integration. It is a monotonically decreasing function of R. The red circles (H = 0.1), 
blue squares (H = 1) and green triangles (H = 10) represent the intermediate membrane 
separations respectively. These symbols represent the same H values in the rest of the 
figures in this Communication. 

For H ^> 1, we can approximate coth(KH) w 1 so that the if -dependence drops 
out. In this limit, the integral can be analytically performed to yield 

d£ 1] (R) = -jp + R- [Hl(i?) - Yl{R)] ' 

where H n (i) are the Struve functions and Y n (t) are the Neumann functions. Since the 
membranes are effectively isolated from each other for large H, (JTTJ coincides with the 
analytical expression obtained for a single membrane geometry [TT]. For R <C 1, (TTTT) 
has a logarithmic behavior, i.e., D^ l \R <C 1) ~ ln(2/i?) — 7 + 0.5 where 7 = 0.5772 • • • 
is the Euler's constant. In this regime, the outer solvent surrounding the membrane is 
unimportant and the membrane behaves effectively as a pure 2D fluid. In the opposite 
R 3> 1 limit, pip decays algebraically D^\r ^> 1) w 2/R. On recasting in dimensional 
form, we see that this limiting expression is independent of the membrane viscosity 77, 
and dependent only on the solvent viscosity r] s . When H <S 1, on the other hand, 
(TTUT) becomes D ( ^ 1 \r/2)/2. The vanishing thickness of the solvent region ii results in 
a rescaling of v~ l by a factor of two, and hence the resultant expression. The solid and 
dashed lines in figure [2] represent the above limiting cases of large and small H limits, 
respectively. It is thus observed that the presence of the second membrane also has 
a finite contribution to the coupling longitudinal diffusion coefficient within the same 
membrane. 
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Figure 4. Dimensionless £>l (R,H) as a function of i? for various values of H. 
The red circles, blue squares and green triangles correspond to H — 0.1, 1 and 10, 
respectively (color online). The dashed line corresponds to the small H limit. 



Following the same argument, the transverse coupling diffusion coefficient between 
two particles within the same membrane is calculated according to 

1 + 2K + coth(KH) 



poo 

D^ 1] (R,H) = 2 dK 
Jo 



X 



Jo(KR) - 



g(K,H) 

Ji(KR) 



(12) 



KR 

In figure El the variation of D^ 1 \r,H) as a function of R for various H is plotted 
using numerical integration. Similar to the longitudinal case, it is also a monotonically 
decreasing function of R. For H ^> 1, (fT2l) can be analytically treated to have the form 



D { T n) (R) = ^-^[Hi(R)-Y 1 (R)] 

+ 7r[U (R)-Y {R)), (13) 

which also coincides with the expression for the transverse coupling diffusion coefficient 
in a single membrane [TT]. In this case, the asymptotic behaviors for small and large R 
are D^ 1] (R < 1) w \n(2/R) - 7 - 0.5 and D^\R » 1) w 2/R 2 , respectively. In the 
opposite limit of H <C 1, f lT2"j) becomes D^ 1 \r/2)/2. As before, the solid and dashed 
lines in figure [3] show the limiting cases of large and small H, respectively. 

Now we proceed to calculate the coupling diffusion coefficients of two particles 
in different membranes. In this case, the longitudinal coupling diffusion coefficient 
expressed in dimensionless units is 

D? ), fl , ff)=2 r dif ^^M. (14) 



g(K,H) KR 

The functional dependence of (R, H) on R for various H is shown in figure HI 
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Figure 5. Dimensionless {R,H) as a function of R for various values of H. 

The red circles, blue squares and green triangles correspond to H — 0.1, 1 and 10, 
respectively (color online). The dashed line corresponds to the small H limit. 



Since cosech^if ) « for H ^> 1, the above integral vanishes for large inter- membrane 
distances as seen by the green triangles (H = 10) in figure SJ This is a reasonable result 
as the membranes are effectively independent of each other. In the limit of H <C 1, (fl4|) 
results in D^ 1 \r/2)/2 which is plotted by the dashed line in figure HI As mentioned 
earlier, this function has an initial logarithmic behavior followed by an asymptotic 1/R- 
decay. 

The transverse coupling diffusion coefficient between two particles in different 
membranes is given by 

cosech(KH) 



poo 

D% 2) (R,H) = 2 dK 
Jo 



g(K,H) 

MKR) - . (15) 

In figure El the variation of D^ 2) (R, H) function of R for various H is plotted. The 
above integral also vanishes when H S> 1 as expected for the decoupled membranes. In 
the opposite limit of H 1, ([15]) results in D^ 1 \r/2)/2 as plotted by the dashed line 
in figure [H 

From figures H] and El it can be observed that for large H, the coupling diffusion 
coefficients vanish. This is due to the exponential decay of the cosech(i^if) term. Up 
to H = 1, the proximity to the second membrane leads to additional contributions to 
the coupling diffusion coefficients. This is the main result of this work. Using typical 
values for the solvent (water) and membrane (lipid bilayer) viscosities, the hydrodynamic 
screening length z/" 1 can be estimated to be of the order of jum. This implies that an 
adjacent membrane within this distance (H < 1) can have a strong bearing on the 
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diffusion dynamics such as for typical stacked supported membrane experiments [9]. 
Qualitatively, the presence of the second membrane can enhance the effective coupling 
diffusion. 

Even though the model presented in this Communication captures the essential 
physics, it looks somewhat simplistic in several respects. We have neglected the finite 
size effect of the membrane inclusions which are known to modify the membrane response 
at small inter-particle distances [TTj 112] . At distances much larger than the inclusion 
size, however, these effects become unimportant. Curvature effects of the membrane is 
significant when the radius of curvature is of the order of the hydrodynamic screening 
length v~ l [TT1 fl8] . The extension of our work to membranes in spherical configuration 
is possible. The out-of-plane fluctuations of the membrane lead to a reduction in the 
diffusion coefficient of proteins in the single membranes [T9j [20] . However, it is known 
that the presence of a substrate or the second membrane would suppress the out-of-plane 
membrane fluctuations [21j[22]. Overall, we expect that fluctuations and the presence 
of a substrate will not qualitatively affect our results. 

To summarize, we have calculated the longitudinal as well as transverse coupling 
diffusion coefficients of particles undergoing Brownian motion in a stacked double 
membrane geometry. We obtained the additional contribution to the coupling diffusion 
coming from the proximity to adjacent membranes. As in biological systems such 
as tissues, the cells are rarely isolated, our results also imply that the proximity to 
neighboring cells can affect the diffusion of objects in cell membranes. 
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